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A predictive model of spray combustion nust incorporate models for the wide 
variety of physical environments in a practical combustor. In regions where 
droplets are closely spaced, combustion resembles a diffusion flame-, where they 
are well-separated, an envelope or wake flame results. The relative velocity 
field between the fuel droplets and oxidizer influences boundary layer develop- 
ment about the droplet, recirculating flow patterns, and droplet shape and 
stability. A model mist encorpass these interacting tenporal and spatial 
effects as well as corcplicated combustor boundaries. The objective of the 
current work is to develop the triangular gridding method for describing the 
individual and collective properties of vaporizing and burning fuel droplets. 

Our approach to this problem has been to modify the basic two dimensional 
Lagrangian model [l] to simulate flows in and about fuel droplets [2]. This 
general-connectivity triangular grid permits accurate representation of 
boundaries as well as material surfaces and interfaces. It also allows variable 
resolution through the insertion of new cells as required to maintain accuracy. 
The finite-difference operators for divergence, curl and gradient are 
constructed to exactly reflect the properties of the continuum operators. The 
construction assures conservation of vorticity and mass and provides a 
determination of the local grid connectivity based on convergence criteria for 
the solution of Poisson’s Equation. Extensions of the model have been made in 
the same spirit: finite-difference operators conform to the continuum limit and 

physical quantities are conserved properly. An important factor to note is that 
development of this method is entirely original work; there is no basic lore to 
fall back on when something goes wrong, and algorithms for a particular type of 
term mist be devised and the best one chosen. 

To date the basic hydrodynamic code incorporates algorithms which allow us 
to include the effects of surface tension and viscosity. These algorithms have 
been incorporated and tested extensively as described below. We have also 
devised algorithms for including the effects of compressibility for subsonic 
flow and for incorporating the effects of thermal conductivity. While these 
latter effects are being tested and incorporated, we will be developing the 
algorithms for vaporization and molecular diffusion. The final steps will allow 
us to describe a burning droplet. 

The test problems performed have included sirailations of inconpressible 
flows about droplets in which the density ratios of droplet to background mate- 
rial have been 2:1, 10:1, and most recently 800:1. The two lower density ratios 
were used to test the various algorithms. The latest tests are aimed at model- 
ling kerosene in air and have incoporated the egressions for surface tension 
and viscosity. 

As an exanple of the way in which development of algorithms has proceeded, 
we describe the latest development in an algorithm for viscosity. Our viscosity 
algorithm originally ewressed the change in the vorticity, at a grid point 
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due to viscosity by d^/dt = where v is the viscosity coefficient. All 

triangle velocities about a vertex contributed equally to the change in vorti- 
city at the central grid point. Although the algorithm produced the correct 
spreading rates for the test case of a shear profile, we found that it only did 
so for fairly regular grid geometries because of the ambiguity in determining 
how the changes in vorticity are translated to velocity changes for different 
grid geometries. For an arbitrary grid, a more complete prescription was neces- 
sary and a method was developed in which V^V is a triangle-centered 
quantity. 

This new triangle -centered algorithm was tested in a calculation of the 
spreading of a shear layer of initially zero thickness. The way in which the 
velocity distribution across this layer evolves and the growth of the width of 
the layer are known quantities which may be coiripared to the results of calcula- 
tions. We found that the calculated layer width agreed exactly with the theory 
and the shear layer velocities were correct over the entire mesh. The 00115 ) 0 - 
nents of velocity perpendicular to the shear layer remained zero indicating that 
the algorithm worked well even for the distorted grid used in the test prob- 
lems. 


In order to test our algorithm for surface tension, we performed calcula- 
tions of droplets which oscillate under the effects of surface tension. The 
results of sinulations could be conpared to the linear theory for small anpli- 
tude oscillations on cylindrical Jets, as first given by Rayleigh. We extended 
this theory to predict droplet oscillation frequencies for a droplet in a back- 
ground gas of finite density. 

In the numerical calculations we studied an n = 2 oscillation for a droplet 
density of 2 g/cm^ and an external fluid density of 1 g/cm^. From the calcula- 
tions shown in Figure 1 , we find that the n\americal oscillation period is 
approximately 1.25xl0“^s whereas the theoretical period is 1.13><10"^s. Most of 
the small discrepancy between the numerical and theoretical results can be 
explained by the finite grid s^'^^cing. However, given Rayleigh’s experience with 
large anplitude oscillations, it is reasonable to expect our conputational 
period to differ somewhat from that given by the linear theory. 

Figures 2 and 3 show other early test calculations. Figure 2 illustrates 
the case for which the density ratio is 2:1 and there is no surface tension or 
viscosity present. We see that a recirculation zone forms early in the calcula- 
tion, compressing the droplet in the direction parallel to the flow. Flow 
within the droplet is initiated by this conpression in a direction normal to the 
external flow. The bulges formed at the top and bottom of the distorted droplet 
are pulled around the recirculation zone by the shear flow which is at a maximum 
at these points. The internal droplet flow is therefore driven by the corrpres- 
sion set up between the front and rear stagnation points and by the high shear 
flow which extends around the top and bottom of the droplet and recirculation 
zone. Ihe interaction of the droplet back onto the external flow occurs pri- 
marily throiigh the enlarged cross-sectional area of the droplet which increases 
the size of the recirculation zone. Eventually the droplet is squeezed into a 
thin layer coating the recirculating zone. The thinned film then shatters into 
several smaller pieces, first at the rear of the droplet and later in the more 
laminar flow toward the front on the droplet. 
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Figure 3 shows the results of a calculation with surface tension for the 
same initial conditions as used in the calculation without surface tension 
(Figure 2). As in the case without surface tension, the internal droplet flow 
is driven by conpression parallel to the external flow and is initially norml 
to the external flow. A recirculation zone is formed in the wake of the 
conpressed droplet and the droplet is deformed as it is swept outward and back- 
ward by the external flow both outside and inside the recirculation zone. 
However, the presence of surface tension provides counteracting forces at 
regions of high curvat\ire. Such forces at the sides and rear of the droplet are 
sufficient to stop the droplet from thinning around the recirculation zone. 

The presentation will summarize the calculations shown above and then 
proceed to describe the more recent calculations at an 800;1 density ratio. 
Introducing this large ratio initially caused several problems in defining 
velocities when cells near the droplet boundary were divided. Ihis problem was 
due to finite resolution effects at the boundary and has been fixed. Current 
simulations with viscosity should give us a realistic picture of flows in the 
droplet itself. 

Finally, we describe the new algorithms to be tested for subsonic conpres- 
sibility and thermal conductivity. A new parametric representation of taut 
splines will also be presented. Ihis extension was necessary to achieve a 
”wiggle-free” spline fit in smooth regions of a droplet interface which were 
near surface discontinuties. 
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Figure 2 
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